#include "AnaDef.h"
#include "XGLUtils.h"
#include "AnaUtils.h"

const Double_t feta = 0.9;

Double_t invcor(const Int_t itype, const Double_t ptmean)
{
  const Double_t dypion = XGLUtils::DeltaY(AnaUtils::Mass(AnaUtils::kPION), ptmean, 0, feta);
  const Double_t dytype = XGLUtils::DeltaY(AnaUtils::Mass(itype), ptmean, 0, feta);

  return dytype/dypion;
}

void copyGr(const Int_t itype, const TGraphAsymmErrors * gfrom, TGraphAsymmErrors * gto, TGraphAsymmErrors * gc)
{
  for(Int_t ipt = 0; ipt<gfrom->GetN(); ipt++){
    const Double_t xfrom = gfrom->GetX()[ipt];
    const Double_t fac = invcor(itype, xfrom);

    const Int_t np = gto->GetN();
    gto->SetPoint(np, xfrom, gfrom->GetY()[ipt]*fac);
    gto->SetPointError(np, gfrom->GetEXlow()[ipt], gfrom->GetEXhigh()[ipt],  gfrom->GetEYlow()[ipt]*fac,  gfrom->GetEYhigh()[ipt]*fac);

    gc->SetPoint(np, xfrom, fac);
  }
}

TGraphAsymmErrors * combGr(const Int_t itype, const TGraphAsymmErrors * gr1, const TGraphAsymmErrors * gr2, TGraphAsymmErrors * & grcor)
{
  TGraphAsymmErrors * gout=new TGraphAsymmErrors; 

  grcor = new TGraphAsymmErrors;

  copyGr(itype, gr1, gout, grcor);
  copyGr(itype, gr2, gout, grcor);

  return gout;
}

void ReadGr(const TString tag, TGraphAsymmErrors * grstat[],  TGraphAsymmErrors * grsys[])
{
  const TString spar[]={"Proton","","Kaon",""};
  TString stmp;
  TH1 * htmp=0x0;
  for(Int_t itype=0; itype<4; itype++){
    grstat[itype] = 0x0;
    grsys[itype]  = 0x0;
    if(spar[itype].Length()==0){
      grstat[itype] = new TGraphAsymmErrors; 
      grsys[itype] = new TGraphAsymmErrors;
    }
    else{
      stmp = Form("%sToPionStat_%s_7000GeV", spar[itype].Data(), tag.Data());
      htmp = (TH1*) gDirectory->Get(stmp);
      if(!htmp){
        printf("no grstat[itype]  itype %d %s\n", itype, stmp.Data()); exit(1);
      }
      grstat[itype] = XGLUtils::HistToGraph("stat",htmp);

      stmp = Form("%sToPionSyst_%s_7000GeV", spar[itype].Data(), tag.Data());
      htmp = (TH1*) gDirectory->Get(stmp);
      if(!htmp){
        printf("no grsys[itype]  itype %d %s\n",itype, stmp.Data()); exit(1);
      }
      grsys[itype] = XGLUtils::HistToGraph("sys",htmp);
    }

    grstat[itype]->SetName(Form("grstat%d", itype));
    grsys[itype]->SetName(Form("grsys%d", itype));
  }
}

int main(int argc, char * argv[])
{
  for(Int_t ii = 0; ii<argc; ii++){
    printf("%d: %s\n", ii, argv[ii]);
  }

  if(argc!=2){
    printf("argc!=2\n"); exit(1);
  }

  TFile::Open(argv[1]);


  //===
  TGraphAsymmErrors * grHIGHstat[4];
  TGraphAsymmErrors * grHIGHsys[4];
  ReadGr("HighpT", grHIGHstat, grHIGHsys);

  TGraphAsymmErrors * grLOWstat[4];
  TGraphAsymmErrors * grLOWsys[4];
  ReadGr("LowpT", grLOWstat, grLOWsys);

  //===

  TList *lout=new TList;

  TGraphAsymmErrors * grstat[4];
  TGraphAsymmErrors * grsys[4];

  TGraphAsymmErrors * grcorstat[4];
  TGraphAsymmErrors * grcorsys[4];

  for(Int_t itype=0; itype<4; itype++){
    grstat[itype]=combGr(itype, grLOWstat[itype], grHIGHstat[itype], grcorstat[itype]);
    lout->Add(grstat[itype]); grstat[itype]->SetName(Form("grstat%d",itype));
    lout->Add(grcorstat[itype]); grcorstat[itype]->SetName(Form("grcorstat%d",itype));

    grsys[itype]=combGr(itype, grLOWsys[itype], grHIGHsys[itype], grcorsys[itype]);
    lout->Add(grsys[itype]); grsys[itype]->SetName(Form("grsys%d",itype));
    lout->Add(grcorsys[itype]); grcorsys[itype]->SetName(Form("grcorsys%d",itype));
  }

  TFile *fout=new TFile(Form("PWGLF_ETA%02.0f.root", feta), "recreate");
  lout->Write();
  fout->Save();
  fout->Close();

  printf("convertPWGLFdonedone!!\n");
  return 0;
}
